library(stress.test.plot.report)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.3.2
# tiltrisk_input_path <- here::here("workspace","tiltrisk","tiltrisk_alpha.csv")
tiltrisk_input_path <- here::here("workspace","tiltrisk","tiltrisk_nrisk.csv")
# Change_NPV_shock
# Change_NPV_with_ecosystem_cost
# Change_NPV_with_ecosystem_social_cost
tiltrisk_trajectories_input_path <- here::here("workspace","tiltrisk","tiltrisk_trajectories.csv")
tiltrisk_df <- readr::read_csv(tiltrisk_input_path) %>%
rename(activity_name=`Activity Name`)
## Rows: 211620 Columns: 35
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (17): companies_id, company_name, country, main_activity, clustered, act...
## dbl (18): Change_NPV_shock, Change_NPV_with_ecosystem_cost, Change_NPV_with_...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
tiltrisk_df <- tiltrisk_df %>%
mutate(
pd_difference = pd_shock - pd_baseline,
crispy_perc_value_change=Change_NPV_with_physical_risk
)
get_tiltrisk_plot_df <- function(tiltrisk_df, group_cols){
tiltrisk_plot_df <- tiltrisk_df %>%
group_by_at(c(group_cols, "term")) %>%
summarise(
crispy_perc_value_change = mean(crispy_perc_value_change),
pd_baseline=mean(pd_baseline),
pd_shock=mean(pd_shock),
.groups="drop"
) %>%
mutate(
pd_difference = pd_shock - pd_baseline
)
return(tiltrisk_plot_df)
}
basic_group_cols <- c("run_id", "baseline_scenario", "shock_scenario", "country")
tiltrisk_plot_df_companies_main_activity <- get_tiltrisk_plot_df(tiltrisk_df, c(basic_group_cols, "main_activity"))
tiltrisk_plot_df_companies_clustered <- get_tiltrisk_plot_df(tiltrisk_df, c(basic_group_cols, "clustered"))
tiltrisk_cluster_activity <- get_tiltrisk_plot_df(tiltrisk_df, c("run_id", "baseline_scenario", "shock_scenario", "country", "main_activity", "clustered"))
tiltrisk_plot_df_companies_main_activity %>%
readr::write_csv(here::here("workspace","tiltrisk","tiltrisk_agg_main_activity.csv"))
tiltrisk_plot_df_companies_main_activity %>%
readr::write_csv(here::here("workspace","tiltrisk","tiltrisk_agg_clustered.csv"))
for (run_id_i in unique(tiltrisk_plot_df_companies_main_activity$run_id)){
plot_data <- tiltrisk_plot_df_companies_main_activity %>%
filter(run_id==run_id_i)
shock_scenario = unique(plot_data$shock_scenario)
for (loc in unique(plot_data$country)){
crispy_npv_change_plot <- pipeline_crispy_npv_change_plot(
plot_data |> dplyr::filter(term==1, country==loc),
x_var = "main_activity"
) +
ggplot2::ggtitle(paste(shock_scenario, '-', loc))
print(crispy_npv_change_plot)
}
}
for (run_id_i in unique(tiltrisk_plot_df_companies_main_activity$run_id)){
plot_data <- tiltrisk_plot_df_companies_main_activity %>%
filter(run_id==run_id_i)
shock_scenario = unique(plot_data$shock_scenario)
for (loc in unique(plot_data$country)){
pd_term_plot <- pipeline_crispy_pd_term_plot(
crispy_data_agg = plot_data %>% filter(country==loc),
facet_var = "main_activity"
) +
ggplot2::ggtitle(paste(shock_scenario, '-', loc))
print(pd_term_plot)
}
}
# Calculate the average ROC per company for the second plot
df_averages <- tiltrisk_df |>
dplyr::group_by(
company_name,
country,
run_id,
main_activity,
shock_scenario,
shock_year
) |>
dplyr::summarise(
average_roc = mean(crispy_perc_value_change),
average_pd_diff = mean(pd_difference)
) %>%
rename(
`Average PD Difference` = average_pd_diff,
`Average NPV % change` = average_roc
)
## `summarise()` has grouped output by 'company_name', 'country', 'run_id',
## 'main_activity', 'shock_scenario'. You can override using the `.groups`
## argument.
# Reshape data to long format
long_data <- df_averages %>%
tidyr::pivot_longer(
cols = c(`Average PD Difference`, `Average NPV % change`),
names_to = "metric",
values_to = "value"
)
# Calculate means and confidence intervals
data_summary <- long_data %>%
group_by(
run_id,
country,
shock_scenario,
shock_year,
main_activity,
metric
) %>%
summarise(
mean = mean(value),
se = sd(value) / sqrt(n()),
ci_upper = mean + qt(0.975, df = n() - 1) * se,
ci_lower = mean - qt(0.975, df = n() - 1) * se
) %>%
ungroup()
## Warning: There were 8 warnings in `summarise()`.
## The first warning was:
## ℹ In argument: `ci_upper = mean + qt(0.975, df = n() - 1) * se`.
## ℹ In group 1: `run_id = "32a7d72b-480d-4de3-8522-6bfb40cfe2a8"`, `country =
## "austria"`, `shock_scenario = "IPR2023_FPS"`, `shock_year = 2030`,
## `main_activity = "agent/ representative"`, `metric = "Average NPV % change"`.
## Caused by warning in `qt()`:
## ! NaNs produced
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 7 remaining warnings.
## `summarise()` has grouped output by 'run_id', 'country', 'shock_scenario',
## 'shock_year', 'main_activity'. You can override using the `.groups` argument.
for (run_id in unique(data_summary$run_id)) {
plot_data <- data_summary %>% dplyr::filter(.data$run_id == .env$run_id)
for (loc in unique(plot_data$country)){
p1 <-
ggplot(
plot_data %>% filter(country==loc),
aes(x = main_activity, y = mean, fill = metric)
) +
geom_bar(stat = "identity", position = position_dodge()) +
geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper),
width = 0.2,
position = position_dodge(0.9)
) +
facet_wrap(~metric, scales = "free_y") +
scale_fill_manual(values = c(
"Average NPV % change" = "#5D9324",
"Average PD Difference" = "#BAB6B5"
)) +
scale_y_continuous(labels = scales::label_percent(accuracy = 1)) +
r2dii.plot::theme_2dii() +
theme(
axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "none"
) +
labs(y = "Average Value", x = "Business Unit") +
ggtitle(paste(
loc,
"-",
plot_data[1, "shock_scenario"] %>% pull()
# , "- shock year :",
# plot_data[1, "shock_year"] %>% pull()
))
print(p1)
}
}
for (run_id_i in unique(tiltrisk_plot_df_companies_main_activity$run_id)){
data_plot <- tiltrisk_df %>%
filter(run_id==run_id_i, term==1)
filtering <- data_plot %>%
distinct(country ,main_activity ,companies_id) %>%
group_by(country, main_activity) %>%
summarise(nrow=n()) %>%
filter(nrow>10) %>%
select(-nrow) %>%
ungroup()
data_plot <- data_plot %>% inner_join(filtering)
shock_scenario <- unique(data_plot$shock_scenario)
density_plot <- make_density_plots(data_plot,
numeric_values = "pd_shock",
density_var = "country",
group_variable="main_activity"
) +
ggplot2::ggtitle(paste0("Distribution of PD at shock for ", shock_scenario, "shock scenario"))
print(density_plot)
}
## `summarise()` has grouped output by 'country'. You can override using the
## `.groups` argument.
## Joining with `by = join_by(country, main_activity)`
## `summarise()` has grouped output by 'country'. You can override using the
## `.groups` argument.
## Joining with `by = join_by(country, main_activity)`
for (run_id_i in unique(tiltrisk_df$run_id)) {
data_plot <- tiltrisk_df |>
dplyr::filter(run_id == run_id_i)
agg_analysis_data <- data_plot |>
# dplyr::filter(.data$net_present_value_difference != 0) |>
dplyr::select(.data$company_name, .data$crispy_perc_value_change, .data$pd_difference) |>
dplyr::group_by(.data$company_name) |>
dplyr::summarise(
crispy_perc_value_change = mean(crispy_perc_value_change),
pd_difference = mean(pd_difference),
.groups = "drop"
)
# Sorting categories based on value1 in descending order
plot_data <- agg_analysis_data |>
# sample_frac(0.1) |>
dplyr::arrange(dplyr::desc(.data$crispy_perc_value_change)) |>
dplyr::mutate(company_name = factor(.data$company_name, levels = .data$company_name)) |>
tidyr::pivot_longer(cols = c("crispy_perc_value_change", "pd_difference"), names_to = "variable", values_to = "value")
# Plotting
p1 <- ggplot(plot_data %>% filter(variable == "crispy_perc_value_change"), aes(x = factor(company_name), y = value, group = variable)) +
geom_step(color = "#5D9324", size = 1) +
scale_y_continuous(
labels = scales::percent_format(accuracy = 1),
breaks = scales::pretty_breaks(n = 5)
) +
geom_hline(yintercept = 0, color = "lightgray", linetype = "dashed", size = 0.5) +
r2dii.plot::theme_2dii() +
theme(
axis.text.x = element_blank(), # element_text(angle = 90, vjust = 0.5),
axis.ticks.x = element_blank(),
axis.title.y = element_text(size = 11),
strip.background = element_blank(),
strip.placement = "outside",
legend.position = "none"
) +
labs(x = NULL, y = NULL) +
guides(fill = NULL) +
ylab("Mean company percent value change")
# Function to create bins every 10 observations
bin_data <- function(data, bin_size) {
data <- data %>%
mutate(bin = (as.numeric(row_number()) - 1) %/% bin_size) %>%
group_by(bin) %>%
summarise(
avg = mean(value),
min = min(value),
max = max(value)
) %>%
ungroup()
return(data)
}
# Bin data every 10 observations
binned_data <- bin_data(plot_data %>% filter(variable == "pd_difference"), round(nrow(plot_data) / 100))
# Create the plot
p2 <- ggplot(binned_data, aes(x = factor(bin), y = avg)) +
geom_col(fill = "#BAB6B5") +
geom_errorbar(aes(ymin = min, ymax = max), width = 0.2) +
scale_y_continuous(
labels = scales::percent_format(accuracy = 1),
breaks = scales::pretty_breaks(n = 5)
) +
r2dii.plot::theme_2dii() +
theme(
axis.text.x = element_blank(),
axis.ticks.x = element_blank(),
axis.title.y = element_text(size = 11),
strip.background = element_blank(),
strip.placement = "outside"
) +
labs(x = NULL, y = NULL) +
guides(fill = NULL) +
ylab("Mean climate Transition-related PD difference")
le_plot <- cowplot::plot_grid(p1, p2, ncol = 1, align = "v")
title <- cowplot::ggdraw() +
cowplot::draw_label(
paste(
data_plot[1, "shock_scenario"] %>% pull(),
" - shock year :",
data_plot[1, "shock_year"] %>% pull()
),
fontface = "bold",
x = 0.5,
hjust = 0.5
) +
theme(
# add margin on the left of the drawing canvas,
# so title is aligned with left edge of first plot
plot.margin = margin(0, 0, 0, 7)
)
le_plot <- cowplot::plot_grid(
title, le_plot,
ncol = 1,
# rel_heights values control vertical title margins
rel_heights = c(0.1, 1)
)
print(le_plot)
}
## Warning: Use of .data in tidyselect expressions was deprecated in tidyselect 1.2.0.
## ℹ Please use `"company_name"` instead of `.data$company_name`
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## Warning: Use of .data in tidyselect expressions was deprecated in tidyselect 1.2.0.
## ℹ Please use `"crispy_perc_value_change"` instead of
## `.data$crispy_perc_value_change`
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## Warning: Use of .data in tidyselect expressions was deprecated in tidyselect 1.2.0.
## ℹ Please use `"pd_difference"` instead of `.data$pd_difference`
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
EU LIFE Project Grant
Scientific Transition Risk Exercises for Stress tests & Scenario Analysis has received funding from the European Union’s Life programme under Grant No. LIFE21-GIC-DE-Stress under the LIFE-2021-SAP-CLIMA funding call.